function Table4

addpath('../Data')
addpath('basic_functions')

var1 = 'RV - IV';
Nfactors = 3;

[Return,Z,IsF,Nt,Chars,CharsNames,~,~,~] = package_data_all_greeks;
Return = Return.ret_daily;
[N,T,L] = size(Z);

% construct W and X matrices that will be used for IPCA
W = NaN(L,L,T);
X = NaN(L,T);
for t = 1:T-1
    W(:,:,t) = squeeze(Z(IsF(:,t+1),t,:))'*squeeze(Z(IsF(:,t+1),t,:)) / Nt(t+1);
    X(:,t+1) = squeeze(Z(IsF(:,t+1),t,:))'*Return(IsF(:,t+1),t+1) / Nt(t+1);
    Return(~IsF(:,t+1),t+1) = NaN;
end
clear t;

% run IPCA
[G, F] = IPCA(X,W,Nt,Nfactors);

% compute expected returns with looking ahead bias
EReturn = NaN(N,T,'single');
Betas = NaN(N,Nfactors,T,'single');
Alpha = NaN(N,T,'single');
for t = 2:T    
    Betas(:,:,t) = squeeze(Z(:,t-1,:))*G;
    EReturn(:,t) = Betas(:,:,t)*F(:,t);
    Alpha(:,t) = Return(:,t) - Betas(:,:,t)*F(:,t);
end

% XXX = Chars;
XXX = Z;
RET = NaN(T,10);
ERET = NaN(T,10);
EE = NaN(T,10,Nfactors);
ALPHA = NaN(T,10);
BETAS = NaN(T,10,Nfactors);

f = [find(strcmp(CharsNames,var1)) find(strcmp(CharsNames,'RV - IV'))];
CHAR = NaN(T,10,length(f));
for t = 2:T
    ind = ~isnan(Return(:,t)) & ~isnan(Z(:,t-1,f(1)));
    if sum(ind) > 0
        rr = Return(~~ind,t);
        er = EReturn(~~ind,t);
        xx = Chars(~~ind,t-1,f);
        x0 = Z(~~ind,t-1,f(1));
        bb = Betas(~~ind,:,t);
        prc1 = [-inf prctile(unique(x0),[10:10:90]) inf];
        for j = 1:10
            n1 = (x0>prc1(j)&x0<=prc1(j+1));
            
            RET(t,j) = mean(rr(n1),'omitnan');
            ERET(t,j) = mean(er(n1),'omitnan');
            ALPHA(t,j) = mean(rr(n1)-er(n1),'omitnan');
            BETAS(t,j,:) = mean(bb(n1,:),'omitnan');
            CHAR(t,j,:) = squeeze(mean(x0(n1,1,:),'omitnan'));
            %CHAR(t,j,:) = squeeze(mean(xx(n1,1,:),'omitnan'));
            EE(t,j,1) = mean(bb(n1,1)*F(1,t),'omitnan');
            EE(t,j,2) = mean(bb(n1,2)*F(2,t),'omitnan');
            EE(t,j,3) = mean(bb(n1,3)*F(3,t),'omitnan');
        end
    end
end

XX = mean(RET,'omitnan');
YY = sqrt(sum(~isnan(RET(:,1)))-1)*mean(RET,'omitnan')./std(RET,'omitnan');
bb = [];
for j = 1:size(XX,2)
    bb = [bb '&' num2str(100*XX(1,j),'%1.2f')];
end
disp(['Ret' bb '\\'])

aa=[];

XX = mean(ERET,'omitnan');
YY = sqrt(sum(~isnan(ERET(:,1)))-1)*mean(ERET,'omitnan')./std(ERET,'omitnan');
bb = [];
for j = 1:size(XX,2)
    bb = [bb '&' num2str_math(100*XX(1,j),'%1.2f')];
end
disp(['E[Ret]' bb '\\'])

XX = mean(ALPHA,'omitnan');
YY = sqrt(sum(~isnan(ALPHA(:,1)))-1)*mean(ALPHA,'omitnan')./std(ALPHA,'omitnan');
bb = [];
for j = 1:size(XX,2)
    bb = [bb '&' num2str_math(100*XX(1,j),'%1.2f')];
end
disp(['$\alpha_{IPCA}$' bb '\\'])

disp('\\')
XX = squeeze(mean(EE,'omitnan'))';
for i = 1:Nfactors
    bb = [];
    for j = 1:size(XX,2)
        bb = [bb '&' num2str_math(100*XX(i,j),'%1.2f')];
    end
    disp(['$\beta_{F' num2str(i) '}$ E[F' num2str(i) ']' bb '\\'])
end
disp('\\')

var = CharsNames(f);
XX = squeeze(mean(CHAR,'omitnan'))';
for i = 2:size(var,1)
    bb = [];
    for j = 1:size(XX,2)
        bb = [bb '&' num2str_math(XX(i,j),'%1.2f')];
    end
    disp([var{i,1} bb '\\'])
end

return